Sort data by age group and month
ds<-ds[order(ds$age_group, ds$date),] #Sort data by age group and monthSet parameters for analysis
Here we need to set a few parameters. We Use the evaluatr.init() function to specify the name of the dataset, the date at which the vaccine is introduced, the date at which we want to begin evaluating the vaccine (typically 1-2 year after vaccine introduction). We also provide some information on the dataset, sch as whether the data are monthly or quarterly (n_seasons), the variable names for the grouping variable, the date variable, the outcome variable, and the denominator variable (if any).
analysis <- evaluatr.init(
country = "Mexico_state", data = ds,
post_period_start = "2008-02-01", #First 'post-intervention' month is Jan 2012
eval_period_start = "2009-02-01", #We ignore first 12 month of data to allow for vaccine ramp up
eval_period_end = "2013-12-01", #The evaluation period lasts 2 years
n_seasons = 12, #This is monthly data, so select 12
year_def = "epi_year", # we are in southern hemisphere, so aggregate results by calendar year (Jan-Dec)
group_name = "age_group", #Strata categry name
date_name = "date", #Date variable name
outcome_name = "J12_18", #Outcome variable name
denom_name = "ach_noj" ,#Denominator variable name
set.burnN=10000,
set.sampleN=5000
)
set.seed(1)Run the main analysis
Save the results in object ‘impact_results’
n.iter<-ncol(impact_results$full$rr_iter)
thin=1
trace1<-impact_results$full$rr_iter[,seq(from=1, to=n.iter, by=thin)] #thin
for(i in 1: nrow(trace1)){
geweke.p<- pnorm(abs(geweke.diag(mcmc(trace1[i,]))$z),lower.tail=FALSE)*2
print(geweke.p)
if(geweke.p>0.05){con.stat<-'Model converged'
}else{
con.stat<-'Not converged'
}
plot(trace1[i,], type='l', main=paste0(con.stat,' ',colnames(trace1)[i]), bty='l', ylim=c(0.2,2))
}
#> var1
#> 0.362741
#> var1
#> 0.2256035
#> var1
#> 0.532991
#> var1
#> 0.004133511
#> var1
#> 0.1015424
#> var1
#> 0.7423961
#> var1
#> 0.004454227
#> var1
#> 0.7486577
#> var1
#> 0.2652687
#> var1
#> 0.8571663
#> var1
#> 0.3705752
#> var1
#> 0.2969788
#> var1
#> 0.6601505
#> var1
#> 0.9622489
#> var1
#> 0.3925018
#> var1
#> 0.04765536
#> var1
#> 0.6839096
#> var1
#> 0.4085838
#> var1
#> 0.3999231
#> var1
#> 0.5864752
#> var1
#> 0.925462
#> var1
#> 0.9224339
#> var1
#> 0.8633828
#> var1
#> 0.5127307
#> var1
#> 0.9975958
#> var1
#> 0.02664334
#> var1
#> 0.182912
#> var1
#> 0.01926025
#> var1
#> 0.3963716
Finally, we can look at the cumulative cases prevented. In this example, there have been 445 cases prevented (95%CrI: 58, 931) from the time of vaccine introduction to the last day month of the study period. This is calculated by takin the difference between the observed and fitted number of cases in each month, and summing them. If atleast 1 control disease is identified from the synthetic controls model, then the result here is drawn from that model, otherwise, it is drawn from the STL+PCA model.
